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pointing  at  bright  objects  that  may  damage  the  equipment.  Additionally,  eigenaxis 
maneuvers  are  the  primary  method  by  which  the  attitude  is  controlled. 

In  this  thesis,  optimal  control  theory  is  applied  to  provide  automated 
maneuver  design  capabilities  to  support  the  LRO  mission.  The  approach  allows 
dynamic  constraints,  as  well  as  other  constraints  such  as  occultation  avoidance, 
to  be  easily  incorporated  into  the  maneuver  design  process.  This  aspect  also 
simplifies  maneuver  checkout  activities.  The  results  of  this  thesis  show  that 
maneuvers  can  be  designed  to  reorient  the  LRO  in  the  presence  of  multiple 
occultation  constraints.  Moreover,  maneuver  times  can  be  reduced  up  to 
90  percent  compared  to  the  conventional  approach.  This  increases  the  potential 
for  efficient  science  data  collection. 
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I.  INTRODUCTION  AND  BACKGROUND 


As  computer  technology  improves,  new  tools  are  continuously  becoming 
available  to  perform  complex  calculations  for  design,  simulation,  and  operation. 
One  field  which  has  largely  remained  unchanged  despite  these  advances  is 
control  theory.  In  control  theory,  handling  additional  complexity  often  requires 
large  matrix  operations  or  iterative  calculations,  something  which  was  not 
feasible  before  low  power,  high  throughput  computers.  In  the  field  of  control 
theory,  particularly  in  control  of  space  systems,  implementation  of  new 
technology  has  been  slow,  and  control  theory  has  tended  towards  those  “things 
that  work”  [1].  Optimal  control  theory  introduces  a  way  to  leverage  these  newly 
available  tools  to  design  and  implement  seemingly  counterintuitive  operations 
which  are  not  viable  using  traditional  methods,  yet  may  improve  performance 
with  respect  to  various  cost  metrics. 

This  objective  of  this  thesis  is  to  explore  the  application  of  optimal  control 
theory  to  the  Lunar  Reconnaissance  Orbiter  (LRO).  Specifically,  the  automation 
of  attitude  maneuver  design  and  checkout  involving  a  real  spacecraft  (S/C)  is 
examined,  using  the  framework  of  optimal  control,  wherein  optimality  is  viewed 
as  a  side  benefit  but  not  necessarily  a  requirement. 

A.  LUNAR  RECONNAISSANCE  ORBITER 

The  Lunar  Reconnaissance  Orbiter  (LRO)  is  a  National  Aeronautics  and 
Space  Administration  (NASA)  mission  launched  on  June  18,  2009.  Its  primary, 
one  year  mission  was  lunar  data  collection  to  include  temperature  mapping,  color 
imaging,  and  ultraviolet  (UV)  albedo  values  [1].  At  the  conclusion  of  the  initial 
mission,  the  LRO’s  mission  was  extended  to  continue  various  science  objectives. 
These  datasets  are  made  publicly  available  in  the  Planetary  Data  System  (PDS) 
[2].  The  collection  of  data  by  the  LRO  paves  the  way  for  future  lunar  mission 
planning  by  NASA  and  others,  including  manned  missions  [1]. 
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The  LRO  payload  consists  of  7  instruments  [1,  2]: 

•  The  cosmic  ray  telescope  for  the  effects  of  radiation  (CRaTER), 
which  characterizes  the  lunar  radiation  environment  and  seeks  to 
measure  radiation  effects  on  a  mock  human  tissue. 

•  The  diviner  lunar  radiometer  experiment  (DLRE),  which  maps 
surface  and  subsurface  temperatures  of  the  moon,  looking  for 
potential  landing  hazards  for  future  missions. 

•  The  Lyman  alpha  mapping  project  (LAMP),  which  seeks  to  discover 
water  ice  in  permanently  shadowed  areas. 

•  The  lunar  exploration  neutron  detector  (LEND),  which  is  used  to 
gather  information  on  the  neutron  radiation  environment  and  to 
analyze  for  evidence  of  water  ice. 

•  The  lunar  orbiter  laser  altimeter  (LOLA),  which  generates  a  high 
resolution,  three  dimensional  surface  map  of  the  moon.  It  also 
identifies  permanently  shadowed  or  illuminated  surfaces. 

•  The  Lunar  Reconnaissance  Orbiter  camera  (LROC),  which  is  a 
black  and  white  surface  imager  with  a  resolution  as  high  as  one 
meter.  A  color  camera  at  100  meter  resolution  is  also  part  of  the 
suite. 

•  The  miniature  radio  frequency  (Mini-RF)  searches  polar  regions  for 
water  ice  using  radar,  and  acts  as  a  demonstration  of  concept  for 
communication  with  Earth  based  ground  stations. 

A  photograph  of  the  LRO,  while  still  in  development,  is  shown  in  Figure  1. 
All  instruments  with  the  exception  of  the  Mini-RF  have  been  installed  and  are 
labeled. 

To  conduct  maneuvers,  the  LRO  uses  four  reaction  wheels  (RW)  arranged 
in  a  tetrahedral  configuration  [3].  By  using  RWs  as  the  primary  source  of  torque, 
fuel  only  needs  to  be  used  for  momentum  dumping  purposes.  Momentum 
accumulation  in  the  RWs  in  a  lunar  orbit  is  primarily  due  to  gravity  gradient 
effects  [4].  The  moon’s  gravitational  field  is  less  uniform  than  the  Earth’s,  and 
thus  has  a  large  effect  on  spacecraft  (S/C)  orbiting  the  moon. 
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Figure  1 .  The  LRO  and  instruments,  after  [2], 


Many  of  LROs  instrumentation  and  sensors  are  optical,  such  as  the  LROC 
or  the  star  trackers  used  for  attitude  determination.  These  instruments  are 
sensitive  to  light  input,  and  could  be  damaged  if  exposed  to  a  bright  source  such 
as  the  sun.  The  mapping  instruments  are  ostensibly  primarily  nadir  pointing,  and 
do  not  usually  need  to  worry  about  light  from  the  sun.  However,  it  is  possible  that 
when  a  maneuver  of  the  S/C  is  required,  sensitive  instruments  may  be  caused  to 
point  at  or  near  the  sun  during  the  maneuver.  This  drives  the  need  to  design 
maneuvers  in  a  way  that  prevents  such  damage  to  instruments,  referred  to  as 
occultation  avoidance  or  obstacle  avoidance.  Since  obstacle  avoidance  for  the 
LRO  is  done  using  a  trial  and  error  type  approach  [5],  maneuver  design  is  not 
autonomous;  that  is,  that  maneuver  design  must  be  done  with  humans  in  the 
loop.  This  motivates  the  desire  to  automate  maneuver  design  [6]. 
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Now  in  an  extended  mission,  the  LRO  continues  to  collect  science  data.  In 
its  first  five  years,  the  LRO  has  collected  “as  much  data  as  all  other  planetary 
missions  combined”  [7],  The  mission  has  been  successful,  and  aims  to  continue 
in  its  legacy  to  better  our  understanding  of  our  closest  neighbor.  By  implementing 
automated  optimal  maneuver  design,  more  of  the  focus  can  shift  from  maneuver 
design  to  the  needs  of  the  science  mission.  This  would  simplify  ground 
operations  and  increase  the  speed  and  volume  of  science  collection. 

B.  OPTIMAL  CONTROL  THEORY 

Control  theory  refers  to  the  branch  of  mathematics  and  engineering  in 
which  the  behavior  of  dynamic  systems  under  certain  input  and  feedback 
conditions  is  analyzed  [8].  Optimal  control  theory,  in  turn,  is  the  mathematics  of 
calculating  an  optimal  strategy  for  the  input  in  order  to  achieve  a  minimized  cost 
[9].  The  definition  of  cost  is  problem  specific,  and  can  be  minimum  time, 
minimum  fuel,  minimum  distance,  etc. 

Pontryagin’s  minimum  principle  provides  the  necessary  conditions  for  an 
optimal  control: 

Pontryagin’s  minimum  principle:  Given  an  optimal  solution 

{.r(-),»(-)T/}  to  Problem  P^,  there  exists  a  costate,  T(.),  and  a 

covector,  v,  that  satisfies  the  Adjoint  Equation,  the  Hamiltonian 

Minimization  Condition,  the  Hamiltonian  Value  Condition,  the 

Hamiltonian  Evolution  Equation,  and  the  Transversality  Condition 

[10,  p.30]. 

The  application  of  Pontryagin’s  minimum  principle  to  optimal  control  theory 
is  discussed  at  length  in  [10],  on  which  this  introduction  is  based. 

To  illustrate  the  application  of  Pontryagin’s  minimum  principle,  a  simple 
example  problem  is  considered.  Consider  a  system  whose  dynamic  model  is  a 
simple  double  integrator,  as  in  (1.1). 


.T  =  V 

v  =  ;/ 


(1.1) 
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The  control  variable  is  u ,  which  is  constrained  in  magnitude  to 
Furthermore,  let  the  problem  be  to  take  the  system  from  an  initial  rest  position  to 
another  rest  position  in  the  fastest  time  possible: 


Minimize; 


X  =  [.r,  v]  G 
u  =  ;/  g[-1,1] 

J[x(-),ii(-)T/]  =  V 


Subject  to: 


.r  =  V 

V  =  // 

(■Yo’Vo)  =  (0,0) 

-!<?/<! 


(1.2) 


In  (1.2),  J  represents  the  cost  function,  which  is  generally  defined  to  be 
the  sum  of  the  endpoint  and  running  costs,  as  seen  in  (1 .3) 

J[x(-),u(-)T^]  =  f(x(V),V)+fV(x(0, 11(0,0^/  (1-3) 

»«o 


The  first  step  in  solving  the  problem  is  to  calculate  the  Hamiltonian.  The 
Hamiltonian  of  (1.2)  is  given  by  (1.4).  In  general,  the  Hamiltonian  is  also  a 
function  of  time.  However,  this  is  not  the  case  in  (1.2)  and  thus  the  time  variable 
is  omitted  from  (1.4)  for  clarity. 

H  (>.,x,u)  =  F(x,u)  +  X^/(x,u) 

^([^x  v],[?/])  =  \v  + 


In  the  Hamiltonian,  the  costate  vector  X  is  introduced.  The  costate  is  a 
vector  associated  with  the  dynamics,  and  is  related  to  the  state  variables  by  the 
adjoint  equation,  given  by 


(1.5) 
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The  Hamiltonian  minimization  condition  is,  literally,  minimization  of  the 
Hamiltonian.  That  is,  the  optimal  control  must  minimize  the  Hamiltonian  function 
over  the  entire  time  horizon.  The  Hamiltonian  minimization  condition  is  given  by 

Min  //([4,/IJ,[.x:,v],[m])  =  4v  +  /1^m 


Subject  to: 


h^<h{u)<h^  (1.6) 


Since  the  control  signal  belongs  in  a  space  bounded  by  inequalities,  h{u) 
the  Lagrangian  of  the  Hamiltonian  must  be  defined  as 


H[]  =  H[-]  +  juh(u) 
H  (ju,  X,  x,u)  =  +  juu 


(1.7) 


The  last  term  includes  a  covector  function  for  the  control  variable.  The 
variable  //  is  a  time  dependent  Karush-Kuhn-Tucker  multiplier  function  which  is 
associated  with  the  constraints  on  the  control  variable  [10].  The  Karush-Kuhn- 
Tucker  condition  states  that  (1.7)  must  be  stationary  with  respect  to  the  control. 
This  gives  the  following; 

^  =  =  0  (1.8) 

ou 


For  every  instant  in  the  time,  the  pair  must  satisfy  the 

complementarity  condition,  which  is  given  by 

<  0,  h{u)  = 

=  0,  h^<h{u)<h^ 

/^i  „  (1-9) 

>0,  h(u)  = 

unrestricted,  =  h’^ 

The  complementarity  condition  applied  to  (1.8)  and  expressed  as  a 
condition  on  the  velocity  costate  is  given  by 


>0, 

u{t)  =  -1 

0, 

-1  <  u(t)  <  1 

(1.10) 

<0, 

u(t)  =  1 
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Using  (1.5),  the  costate  histories  are  given  in  (1.11),  where  a  and  b  are 
arbitrary  constants  of  integration. 


=  0^  =  a 

^  =  -at  -  b 


(1.11) 


From  (1.11),  it  is  evident  that  the  position  and  velocity  costates  are 
constant  and  linear,  respectively,  for  the  example  problem.  Assuming  A^  has  a 

non-zero  value  over  some  non-zero  time  interval,  (1.10)  provides  the  basis  for  a 
switching  function.  That  is,  the  control  signal  will  always  be  positive  or  negative 
one,  depending  on  the  sign  of  T, .  The  result  is  summarized  as 


u(t)  =  sgn(at  +  b)  (1.12) 

If  there  are  not  enough  boundary  conditions  specified  as  part  of  problem 

,  the  transversality  condition  may  be  used  to  solve  for  the  additional  conditions 

necessary  to  find  a  unique  solution  to  the  problem.  The  transversality  condition  is 
given  by 

=  (1.13) 

dx, 

where 

E(v,x,  t^.  ):=E{Xf,t^)  +  v^e(x^  ,t^.)  (1.14) 

The  first  term  on  the  right  side  of  (1.14)  is  the  endpoint  cost  from  (1.3). 
The  elements  of  vector  v  are  the  Lagrange  multipliers  which  are  associated  with 
the  constraints  e(x^,r^)  =  0.  The  introduction  of  additional  unknowns  in  the  form 

of  Lagrange  multipliers  allows  any  missing  boundary  conditions  to  be  solved  for  a 
unique  solution  [10].  For  the  example  problem  given  here,  application  of  the 
transversality  conditions  yield 
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(1.15) 


E{v,x,tjr):=  tj  +Vj(x^  -x'^)  +  V2(vy  -v^) 

Kitf)=Vy 

\(tf)  =  y2 

Since  the  example  problem  had  two  endpoint  conditions  specified  for  each 
state  variable,  the  transversality  condition  does  not  provide  any  additional 
boundary  conditions.  However,  since  the  final  time  is  not  fixed,  the  Hamiltonian 
value  condition  must  also  be  examined: 


dtf 


(1.16) 


From  (1.16),  the  value  of  the  Hamiltonian  at  the  final  time  must  be 
negative  one.  Finally,  the  Hamiltonian  evolution  condition  must  be  evaluated. 
This  is  simply  the  time  derivative  of  the  minimized  Hamiltonian,  and  is  given  as 


dH 

dt 


(1.17) 


The  evaluation  is  not  intuitively  apparent  due  to  the  signum  function.  For 
the  example  problem,  the  solution  is  shown  in  (1.18). 

V  =  sgn{at  +  b) 
at  +  b\ 


V  = 


-  +  c 


4  =  a 
=  -at-b 

d{X^v  +  (sgn(at  +  b)) 


H  =■ 


dt 


(1.18) 


d{a 


at  +  b\ 


+  c 


+  {-at  -  b){sgn{at  +  b)) 


dt 


d{at  +  b  +ac-  at  +  b\) 


dt 


=  0 
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Despite  the  presence  of  time  in  the  equation,  the  Hamiltonian  does  not 
actually  explicitly  depend  on  it.  The  fact  that  the  Hamiltonian  is  constant  at 
negative  one  is  true  for  all  minimum  time  problems. 

While  it  is  possible  to  develop  and  apply  the  necessary  conditions  to  solve 
any  optimal  control  problem,  analysis  quickly  becomes  difficult  for  real  world 
problems  where  additional  constraints  and  dynamic  relationships  must  generally 
be  introduced.  Thus,  numerical  methods  and  computers  are  employed.  One  such 
software  used  to  solve  optimal  control  problems  is  DIDO.  DIDO  is  a  MATLAB 
program  which  solves  any  optimal  control  problem,  within  the  bounds  of  proper 
formulation  and  computing  power  [11].  It  allows  complete  specification  of 
dynamics,  boundary,  and  time  variant  conditions,  leading  to  the  ability  to 
construct  problems  incrementally  without  having  to  come  up  with  all  new  code  for 
the  problem. 
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11.  LUNAR  RECONNAISSANCE  ORBITER  DYNAMIC  MODEL 


In  order  to  apply  optimal  control  theory  to  the  LRO,  the  S/C  dynamics 
need  to  be  appropriately  modeled.  The  dynamics  of  a  S/C  are  affected  by 
properties  such  as  geometry,  motion  of  parts,  thermal  expansion,  external 
torques,  change  in  fuel  level  and  distribution,  etc.  Modeling  every  parameter  of 
such  a  complex  system  quickly  makes  the  problem  too  large  to  solve.  Thus,  the 
important  parameters  to  be  modeled  are  selected  and  certain  assumptions  made 
to  maintain  the  scope  of  the  problem  to  a  tractable  size  while  appropriately 
maintaining  fidelity. 

The  attitude  control  problem  is  defined  relative  to  an  inertial  reference 
frame.  This  ensures  conservation  of  angular  momentum.  The  most  important 
simplification  made  in  this  regard  is  that  orbital  velocities  of  the  LRO  around  the 
Moon,  of  the  Moon  around  the  Earth,  the  Earth  around  the  sun  are  disregarded. 
That  is,  the  LRO  will  be  modeled  at  a  slice  in  time  in  its  orbit.  This  allows  the 
body  frame  to  be  treated  as  an  inertial  frame,  and  greatly  simplifies  the 
mathematics.  Additionally,  it  allows  for  celestial  objects  to  be  modeled  as  being 
stationary  in  the  inertial  frame,  simplifying  the  model  further. 

Physical  properties  and  operational  limits  of  the  LRO  were  obtained  via 
personal  correspondence  with  the  LRO  program  office  [3],  as  well  as  open 
literature  such  as  [5].  They  are  listed  in  Table  1.  The  inertia  tensor  is  an  estimate 
of  the  LRO  at  EOL.  The  reaction  wheel  capacities  listed  are  at  their  operational 
limits,  and  represent  a  value  25  percent  lower  than  design  capacity.  The  power, 
omni  antenna,  and  star  tracker  constraints  are  listed  for  completeness,  but  not  all 
are  implemented  as  constraints  for  maneuver  design. 
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Table  1.  LRO  properties  and  constraints,  after  [3]  and  [5]. 


Inertia  tensor  (7 ) 

“  591.47  -20.718  21.459“ 

-20.718  849.93  -28.031 

21.459  -28.031  922.51 

kg  -m^ 

RW  alignment  (z ) 

“0.57357  0.57357  0.57357  0.57357  “ 

0.57923  -0.57923  -0.57923  0.57923 

0.57923  0.57923  -0.57923  -0.57923 

Angular  distance  of  S/C  +z 

axis  from  sun 

63  deg  minimum 

Maximum  angular  velocity  of 

S/C(co_) 

0.1  deg/sec  per  axis 

Maximum  angular  momentum 

inRW(/._) 

60  N-m-s  each 

Maximum  RW  control  torque 

( ^  ) 

0.16  N-m  each 

Power  constraint 

-y  axis  preferred  toward  sun 

Omni  antenna  constraint 

Earth  placed  in  preferred  quadrant  of  S/C  x-y 

plane 

Star  tracker  constraint 

Field  of  view  of  star  trackers  25  deg  from  sun, 

15  deg  from  Moon,  and  15  deg  from  Earth 

The  key  tenet  of  this  thesis  lies  in  the  pointing  constraints.  The  aim  is  to 
perform  simulations  that  involve  avoidance  maneuvers,  so  that  an  approach  can 
be  developed  for  their  automation.  Therefore,  a  maneuver  scenario  must  be 
developed  which  would  violate  a  pointing  constraint  if  no  action  is  taken.  To  do 
this,  first,  one  of  the  pointing  constraints  is  defined  in  the  problem.  As  per  Table 
1 ,  the  sensitive  instruments  located  on  the  +z  face  of  the  S/C  must  stand  off  from 
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the  sun  by  63  degrees.  Therefore,  the  problem  is  defined  such  that  the  S/C  +z 
axis  must  remain  63  degrees  or  more  away  from  the  direction  of  the  sun. 

The  spacecraft  dynamics  must  also  be  defined.  Quaternions  are  chosen  to 
represent  the  S/C  attitude.  Since  the  LRO  uses  RWs  as  its  only  means  of  attitude 
control,  the  combination  of  (2.1)  and  (2.2)  is  used  as  the  dynamic  model  [12]. 
The  controls  are  defined  to  be  the  RW  torques,  as  seen  in  (2.3).  This  assumes  a 
rigid  body  S/C  and  that  relative  motion  of  the  RWs  is  negligible. 

4i  —  —  (K)i<?4  ~  C02<?3  +  C03<72  ) 

Ql  ~  T  ~  ®3<?i ) 

\  (2.1) 

<^3  —  —  (— C0j^2  +  C02<5'i  +  ®3<?4  ) 

Qa  ~  ®3‘?3) 

&)  =  -j\a)x(Ja)  +  Zh)  +  Zh)  (2.2) 

h  =  u  (2.3) 

Next,  the  boundary  conditions  are  defined.  The  required  specifications 
consist  of  an  initial  and  end  attitude  in  quaternions,  initial  and  end  S/C  angular 
rates,  and  an  initial  RW  momentum  vector.  It  is  not  necessary  to  define  a  final 
RW  momentum  vector,  as  it  is  implicitly  specified  via  conservation  of  angular 
momentum  in  the  inertial  frame.  The  initial  attitude  is  defined  to  allow  the  body 
frame  to  coincide  with  the  inertial  frame.  This  allows  the  reconstruction  of  the 
actual  maneuvers  used  on  the  LRO  to  be  done  more  easily,  as  discussed  in 
Chapter  III. 

Lastly,  a  procedure  for  implementing  the  pointing  constraint  needs  to  be 
defined.  To  do  this,  a  path  constraint  can  be  introduced,  which  relies  on  the 
direction  of  the  sun  in  the  inertial  frame.  Since  the  direction  of  the  sensitive 
instruments  is  defined  in  the  body  frame,  the  instrument  vector  must  be 
converted  to  a  vector  in  the  inertial  frame  by  using  a  DCM  derived  from 
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quaternions  as  in  (2.4),  where  is  a  unit  vector  expressed  in  the  body  frame, 
is  the  same  vector  expressed  in  the  inertial  frame,  and  the  quaternions 
describe  the  relative  orientation  of  the  body  frame  with  respect  to  the  inertial 
frame  [12],  This  is  a  coordinate  transformation,  wherein  a  coordinate  system  is 
rotated  about  an  axis  of  rotation  and  angle  described  by  a  quaternion.  Since  the 
S/C  +z  vector  is  [0  0  1]^  in  the  body  frame,  the  S/C  +z  vector  expressed  in  the 
inertial  frame  is  given  by  the  last  column  of  (2.4). 


'^N 


e 


Q\  Qi  Qt,  "*"^4  ^3,^4) 

2(q,q2  -  q^q,  )  -q^  +  q\  -  ql  +  ql 

2(^A+<?2'?4) 


2(?2?3+?1?4) 

-ql-ql+ql+ql 


e 


(2.4) 


By  utilizing  the  definition  of  the  dot  product,  the  angle  between  two  unit 
vectors  S  and  z  can  be  calculated  as  (2.5).  A  graphical  representation  of  this 
relationship  is  shown  in  Figure  2. 

6'  =  cos‘(5.z)  (2.5) 


inertial  frame 
A  / 

zJ 

direction  of  Sun 
(expressed  in  inertial  frame) 


S/C  body  axes 
(expressed  in  inertial  frame) 

Figure  2.  Illustration  of  sun  avoidance  angle  calculation. 
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III.  RECONSTRUCTION  OF  EXISTING  SLEW  ALGORITHM 


To  establish  a  baseline  with  which  optimal  control  results  shall  be 
compared,  the  actual  types  of  maneuvers  used  on  the  LRO  must  be  recreated.  In 
general,  the  LRO  utilizes  the  eigenaxis  slew  to  perform  point  to  point  slews  [5]. 
An  eigenaxis  maneuver  is  one  in  which  the  axis  of  rotation  is  held  constant.  The 
advantage  of  an  eigenaxis  maneuver  is  the  fact  that  the  angle  of  rotation 
between  any  two  attitudes  is  minimized,  but  this  does  not  necessarily  translate  to 
a  minimum  fuel  or  minimum  time  maneuver  [13]-[15].  Hereafter,  an  eigenaxis 
slew  that  takes  the  S/C  from  a  one  defined  orientation  to  another  will  be  referred 
to  as  a  direct  slew  [5]. 

The  maneuver  of  the  LRO  increases  in  complexity  when  a  pointing 
constraint  is  introduced.  For  example,  if  a  direct  slew  causes  sensitive 
instruments  to  point  at  the  sun  at  any  instant,  then  the  maneuver  cannot  be 
performed.  A  trial  and  error  approach  is  used  to  determine  the  appropriate 
maneuver,  as  described  in  [5]  and  explained  below. 

Consistent  with  [5],  the  term  minimum  slew  will  be  used  to  designate  a 
maneuver  such  that  a  specific  vector  in  the  S/C  body  frame,  such  as  the  vector 
representing  the  direction  that  instruments  may  be  pointing,  is  slewed  point  to 
point.  This  is  also  an  eigenaxis  maneuver,  but  it  differs  from  a  direct  slew  in  that 
only  one  vector  of  the  S/C  is  examined  as  opposed  to  the  entire  coordinate 
frame.  A  minimum  slew  allows,  for  example,  changing  the  direction  that 
instruments  are  pointing  without  regard  to  the  specific  orientation  of  the  entire 
S/C. 

In  all  approaches  for  avoiding  violation  of  a  pointing  constraint,  the  slew  is 
broken  up  into  two  segments.  Let  the  direction  that  an  axis  of  concern  is  pointing 
at  the  beginning  of  a  maneuver  be  called  “Point  A”  and  at  the  end,  “Point  B.”  The 
path  between  Point  A  and  Point  B  represents  the  direction  that  the  axis  of 
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interest  points  during  the  maneuver.  In  an  eigenaxis  maneuver,  this  path  will  be  a 
portion  of  a  great  circle.  Let  the  midpoint  of  this  arc  be  denoted  as  “Point  C.” 

The  first  attempt  to  break  up  the  maneuver  utilizes  a  minimum  slew  of  the 
axis  of  interest  from  Point  A  to  Point  C,  and  then  a  direct  slew  to  take  the  S/C  to 
its  desired  final  orientation.  If  this  maneuver  still  causes  a  pointing  constraint 
violation,  then  the  order  is  swapped.  That  is,  the  direct  slew  is  performed  first. 

If  the  first  method  of  mitigation  fails,  the  next  step  is  to  move  Point  C  away 
from  the  center  of  the  obstacle  by  30  degrees.  The  combinations  of  minimum  and 
direct  slews  are  then  investigated.  If  the  constraint  is  still  violated.  Point  C  is  now 
moved  60  degrees  away  from  the  center  of  the  obstacle,  and  the  process 
repeated.  The  two  stage  maneuvers  will  hereafter  collectively  be  called  “dogleg” 
maneuvers. 

To  illustrate  the  application  of  this  idea,  the  initial  and  final  attitudes  as  well 
as  the  direction  of  the  sun  are  chosen  such  that  a  direct  slew  will  cause  the  +z 
axis  of  the  spacecraft  to  point  toward  the  sun.  This  will  violate  the  standoff 
constraint  for  sensitive  instruments  (See  Table  1).  The  specific  attitude  values  for 
this  example  maneuver  are  detailed  in  Table  2.  When  the  aforementioned 
protocol  is  followed  to  mitigate  this  violation,  the  first  mitigation  method  to  yield  a 
satisfactory  solution  is  to  move  Point  C  by  30  degrees  from  its  original  location, 
and  then  to  perform  a  minimum  slew  followed  by  a  direct  slew.  The  original  and 
modified  maneuvers  are  shown  in  Figure  3. 


Table  2.  Boundary  conditions  for  example  sun  avoidance  maneuver. 


[0,0,0,  If 

[0.400, -0.800,0.4279,0.130f 

(Dq  =  CO^ 

[0,0,0]^  deg/ sec 

16 


Current  LRO  Maneuvers 


-S 


Figure  3.  Direct  slew  and  dogleg  maneuver. 


Figure  3  shows  the  motion  of  the  +z  axis  of  the  spacecraft,  as  well  as  a 
cone  depicting  a  63  degree  half  angle  from  the  direction  of  the  sun.  Thus,  if  the 
path  of  the  +z  axis  goes  inside  the  projection  of  the  red  cone,  representing  the 
standoff  range,  there  is  a  constraint  violation.  The  dogleg  maneuver  is  clearly 
able  to  avoid  the  pointing  violation  that  occurs  with  a  simple  direct  slew.  The 
numerical  values  of  the  angle  between  the  S/C  +z  axis  and  the  sun  are  shown  in 
Figure  4. 


10  20  30  40  50  60  70  80  90  100 

Percentage  of  slew  complete 

Figure  4.  Angular  proximity  of  +z  axis  to  center  of  sun  avoidance  cone. 
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A. 


BASELINE  MANEUVER  PERFORMANCE 


In  order  to  quantify  the  performance  difference  between  the  existing  and 
optimal  control  maneuvers,  the  time  it  takes  to  complete  the  existing  maneuver 
will  be  estimated.  To  simplify  the  calculation,  angular  acceleration  is  assumed  to 
be  infinite.  This  allows  for  the  calculation  to  be  reduced  to  an  angular  distance 
and  velocity. 

In  a  quaternion,  the  first  three  elements  encode  a  unit  vector  e,  which  is 
the  axis  of  rotation,  while  the  last  element  represents  the  angle  of  rotation  of  a 
frame  about  this  axis,  as  shown  in  (3.1).  The  axis  of  rotation  is  broken  down  into 
components  expressed  in  the  inertial  frame. 


(1  = 


.  .  9 
e,  sin(-) 

sin(-) 

e,  sin(-) 

cos(^) 


(3.1) 


Let  q  and  q  represent  consecutive,  separate  coordinate  rotations.  Then 
the  equivalent  quaternion  q  can  be  found  using  equation  (3.2)  [12].  The  equation 
has  been  arranged  such  that  the  quaternion  for  the  second  coordinate 
transformation  is  isolated  for  path  recreation. 


r 

<N 

1 

1 _ 

-1 

ft 

-q[  q\  q[  q'^ 

qi 

Qi  -Q[  ^4  ^3 

ft 

_-9.[  -Qi  -ft  ft^ 

_ft_ 

(3.2) 


Since  the  initial  conditions  of  the  problem  dictate  that  the  S/C  body  frame 
be  coincident  with  the  Inertial  frame,  the  relationship  given  in  (3.1)  allows  a 
simple  calculation  of  the  angle  of  rotation  for  the  direct  slew.  Furthermore,  the 
angular  velocity  of  the  S/C  about  each  axis  must  be  proportional  to  the 
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components  of  the  axis  of  rotation  expressed  in  the  body  frame,  as  in  (3.3).  An 
easy  way  to  visualize  this  is  to  simply  assume  a  rotation  about  the  inertial  +X 
axis.  Since  the  body  frame  and  inertial  frame  are  originally  coincident,  it  follows 
that  the  S/C  must  rotate  about  its  +x  axis  for  this  maneuver.  By  setting  the  axis 
with  the  largest  sectorial  component  to  0.1  degrees/sec,  the  maximum  angular 
velocity  for  the  slew  can  be  calculated. 


where 


0.1  deg/  5 

max{ej  ,  ^2  5  4|} 


(3.3) 


For  the  dogleg  maneuver,  it  is  necessary  to  find  the  new  location  of 
Point  C.  First,  the  original  Point  C  is  found  by  finding  the  direction  of  the  S/C  +z 
axis  at  the  midpoint  of  the  direct  slew.  This  is  done  by  using  (3.1),  using  half  of 
the  original  angle  for  the  fourth  component  and  normalizing  to  ensure  that 

4 

=1.  This  quaternion  is  then  substituted  into  (2.4)  to  locate  the  S/C  +z  axis  at 

i=\ 

the  midpoint  of  rotation,  which  by  definition  is  Point  C. 

To  rotate  this  point  30  degrees  away  from  the  center  of  the  obstacle,  the 
vector  corresponding  to  Point  C  is  rotated  about  an  axis  perpendicular  to  the 
center  of  the  obstacle  and  Point  C.  Let  C  represent  the  original  direction  of 
Point  C  in  the  inertial  frame,  S  the  direction  of  the  sun,  and  e^the  axis  about 

which  Point  C  will  be  rotated.  Then,  the  axis  of  this  rotation  can  be  found  by 
using  (3.4).  By  again  applying  the  coordinate  transformation  of  equation  (2.4) 
using  30  degrees  as  the  angle  of  rotation,  the  unit  vector  pointing  at  the  new 

Point  C,  heretofore  referred  to  as  C  is  located.  Figure  5  shows  an  illustration  of 
this  process. 
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(3.4) 


5xC 

SxC 


Figure  5.  Process  of  locating  C . 

Since  the  first  leg  is  a  minimum  slew,  the  axis  of  rotation  of  the  S/C  body 
frame  is  found  by  taking  the  cross  product  of  the  vectors  corresponding  to  the 
initial  and  final  directions  of  the  S/C  Z-axis,  expressed  in  the  inertial  frame.  For 
the  second  leg,  the  relationship  in  (3.2)  is  used  to  calculate  the  quaternion 
representing  only  the  second  leg  of  the  maneuver.  The  axis  of  rotation  for  the 
second  maneuver  is  expressed  in  the  S/C  body  frame  using  equation  (2.4),  after 
which  equation  (3.3)  is  used  to  determine  the  S/C  angular  velocities. 

The  angular  lengths  and  maneuver  times  of  the  baseline  direct  slew  and 
dogleg  maneuvers  are  detailed  in  Table  3.  For  further  analysis  and  comparison, 
the  total  time  for  the  dogleg  maneuver  is  used,  which  is  1,975  seconds.  Note  that 
the  dogleg  maneuver  takes  48  percent  longer  than  the  direct  slew,  which  could 
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mean  that  an  opportunity  for  science  collection  could  be  lost.  This  motivates  the 
desire  to  perform  a  sun  avoidance  slew  in  the  least  amount  of  time. 


Table  3.  Estimates  of  status  quo  maneuver  times. 


Maneuver 

Angular  Length  (deg) 

Time  (sec) 

Direct  slew 

165 

1,331 

Dogleg  leg  1  (minimum) 

74 

764 

Dogleg  leg  2  (direct) 

128 

1,211 

Dogleg  (total) 

202 

1,975 
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IV.  AUTOMATED  MANEUVER  DESIGN  AND  CHECKOUT 


The  dynamic  model  constructed  in  Chapter  III  as  well  as  the  constraints 
introduced  are  now  used  to  develop  an  optimal  control  problem  for  automated 
maneuver  design.  This  chapter  introduces  the  optimal  control  model  and 
adjustments  that  were  made  to  ensure  proper  formulation  was  present  to  obtain 
flight  ready  solutions.  Automation  is  used  in  the  sense  that  an  avoidance 
maneuver  can  be  designed  without  resorting  to  the  heuristic  steps  described  in 
Chapter  III.  Since  the  maneuver  design  is  obtained  as  the  solution  to  an  optimal 
control  problem,  human  intervention  is  not  required  to  generate  the  slew. 
Moreover,  the  maneuver  will  also  automatically  satisfy  the  vehicle’s  dynamic 
constraints,  a  feature  that  can  be  leveraged  to  help  facilitate  the  pre-flight 
maneuver  checkout. 

A.  OPTIMAL  CONTROL  PROBLEM 

From  the  conditions  described  in  Chapter  II,  the  optimal  control  problem 
statement  becomes: 

X  =  [q,CL),h]  G 
u  G  St"* 

Minimize:  J[x(-),u(-),tf]  =  t^ 
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Subject  to: 


^^1  =  2(“i^4-®2^/3+“3^2) 
^2=^(<»l?3+“2^4-«'3<?l) 

(h  =  ^(-®1^2  +  “2^1  +  ®394) 
94  =  -  “2^2  +  “3^3) 


(Si--J  * (r<> X {J CO + Z//)  +  Z//) 
h  =  u 


«/)=[<?/, ®/] 

|/;,  I  <  60 -vi-s 
|4|<0.16A^-w 
\(o.\  <  O.ldeg/ 5 
1^1  >  63  deg 


(4.1) 


where  the  sun  avoidance  angle  is  given  as: 


6  -  cos  Hs»f) 


In  (4.1),  the  two  vector  quantities  represent  the  direction  of  the  sun  (i) 
and  the  direction  of  the  S/C  +z  axis  (f ),  both  expressed  in  the  inertial  frame. 


1.  Necessary  Conditions 

The  Hamiltonian  for  problem  (4.1)  is  given  by  (4.2).  The  variable  Z  refers 
to  the  costate  variable,  with  the  corresponding  state  variable  listed  as  a  subscript. 
Due  to  the  cross  product  terms  present  in  the  dynamic  equations,  there  are 
myriad  terms  in  the  Hamiltonian,  and  it  is  almost  impossible  to  determine  by 
inspection  the  effects  that  the  different  variables  may  have  on  the  Hamiltonian. 
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H  = 

+^[-4,  ((y„y„  -4,4,,)(#4)-(y„y,.  -y,y„)(#3)+(y„y,  -y„y„)(#2)) 

+2^  ((4„4„  -4,/„)(#4)-(7„7,,  -J,.Jjm)  +  I.J,J,,  -■/„/„)(#2)) 

-/I  ((/  /  -/  /  )(#4)-(J  J  -J  J  )(#3)  +  (J  J  -J  J  )(#2))1 

+  — (<^2^3  “‘?3®2  ‘?4®1  )  “*"  A2  (‘?3®1  ‘?4®2  ) 

“*"^3  (‘?1®2  “  ‘?2®1  ‘?4®3  )  “  A4  ^  ‘?2®2  <?3®3 


(4.2) 


where 


#1  =  7  7  7  -7  7  7  -7  7  7  +7  7  7  +7  7  7 

AX  yy  zz  XA  jz  zy  xy  yx  zz  Ay  yz  ^A  xz  yx  zy 


^  2  =  5)  —  (fl  6)  +  Zj  jMj  Z12W2  ^13^3  ^14^4 

^  3  =  7)  —  5)  +  Z2jM|  +  ^22^2  +  Z23M3  +  Z24W4 

H-  A  —  COy  (Jt  6)  —  (^2  2)  +  +  Z32W2  ^33^3  ^34^4 

#5  —  7j^®j  +  J  ^y<X>2  ■*■  ■^ZZ®3  A  A  A2A  A3A  A4A 

#  6  =  J  yjj)^  +  J  yyCO^  +  ^  y^^^  +  ZjjZZ]  +  ^22^2  +  ^3^  ^4  ^ 

#7  =  7^0j  +  7^®2  +  Az®3  +  a  a  '*■  A2A  A3A  A4A 


-7  7  7 

xz  >7 


ZX 


Hamiltonian  minimization  yields  switching  functions  for  all  four  controls, 
since  these  appear  linearly.  When  the  first  switching  function  is  positive,  the  Mj  is 

the  maximum  negative  value,  and  vice  versa.  Similar  results  follow  for  the  other 
controls.  The  adjoint  equations  are; 
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\=^(-\®1-4®2-\®3) 

^«,=^K-?2-\-?3-^,.-?4  +  \-?i)  +  ^[^«,(#16#5  +  #13#3-#10#17) 
-/l^(#15#5  +  #12#3-#9#17)  +  /l^(#14#5  +  #ll#3-#8#17)] 
2^=iKft-2„9,-4«4+2,.«j)  +  ^[A^(#16#7  +  #10#2-#13#18) 
+4^(#15#7  +  #9#2-#12#18)  +  /l^(#14#7  +  #8#2-#ll#18)" 
2,=i(4«,-\,j-2,_,,-2,^,,)  +  i[2^(#13#6  +  #10#4-#16#19) 
-4^(#12#6  +  #9#4-#15#19)  +  /l^(#ll#6  +  #8#4-#14#19)’ 
=^[4^(#15#31  +  #12#27  +  #9#23)-/l^(#16#31  +  #13#27  +  #10#23) 
-/l^(#14#31  +  #ll#27  +  #8#23)’ 

^'-2  =^[^«^(#15#30  +  #12#26  +  #9#22)-/1^(#16#30  +  #13#26  +  #10#22) 
-/l^(#14#30  +  #ll#26  +  #8#22)' 

A,  =—[4  (#15#29  +  #12#25  +  #9#21)-A  (#16#29  +  #13#25  +  #10#21) 

3  ^IL  2  “h 

-A^(#14#29  +  #ll#25  +  #8#21)’ 

^  =^[^«^(#15#28  +  #12#24  +  #9#20)-A^(#16#28  +  #13#24  +  #10#20) 
-A^(#14#28  +  #ll#24  +  #8#20)’ 
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where 


#\  =  J  J  J  -J  J  J  -J  J  J  +J  J  J  +J  J  J  -J  J  J 

XX  yy  zz  XX  yz  zy  xy  yx  zz  xy  yz  zx  xz  yx  zy  xz  zx 


#2  —  J ~  J ^ zz^i  -^31^  -^32^2  +  Z^Jl^ 

#3  =  2J —  J^COy  +  J ^(02  +  J zz^i  -^31^  -^32^2 
#4  —  J y/O-^  +  J yyC02  zy^^  ^  zz^2  -^22^  -^23^  ^24^4 

#5  —  2J y^CO^  —  J ^0)2  +  J yy(^2  ^ yz^2  -^21^  -^22^  -^23^  -^24^4 

#6  =  +  J xy(^2  "*"  ^■^xz®3  “  "*"  ^1 A  "*"  ^2^2  "*"  ^3^  "*"  A4A 

#7  —  J  xy^l  +  Az®3  ~  yy^l  ~*~  A  A  "*"  ^12^2  "*"  ^3^  "*"  A4A 


zz 

#9 

7 

-7 

7 

>x 

ZZ 

zx 

#10 

=  7 

7 

-7 

7 

zy 

T> 

'  zx 

#11 

-7 

7 

-7 

7 

xy 

ZZ 

xz 

#12 

=  7 

7 

-7 

7 

xx 

ZZ 

xz 

zx 

#13 

=  7 

7 

-7 

7 

xx 

•TV 

zx 

#14 

=  7 

7 

-7 

7 

xz 

TT 

#15 

=  7 

7 

-7 

7 

XX 

xz 

yx 

#16 

A. 

^Jy. 

#\S  =  j^(02-j^yG), 
#\9  =  J^^co2-JyzCO, 


#  20  =  Z24CO2  —  Z^^co2 

#21  =  Z23®3  —  ^33^2 

#22  =  Z  22(02  —Z  22(02 
#23  =  Z2J®3  —  Z3J®2 
#24  =  Zj4(»3  -Z34OJ 
#25  =  Zj3®3  —  Z33®J 

#  26  =  Zj2®3  —  Z32®J 

#27  =  Zjj03  -Z3j(»j 

#  28  =  Zj4®2  “  -224®! 

#  29  =  Zj3®2  “  Z22(Oi 
#30  =  Zj202  “•^22®! 
#31  =  Zjj®2  “Z2J®[ 


The  adjoint  equations  of  (4.4)  can  be  used  to  verify  the  proper  trajectories 
for  the  costates.  Evaluation  of  the  transversality  condition  yields  some 
information  that  will  augment  the  information  in  (4.4)  and  is  useful  for  checking 
optimality.  Specifically,  (r^)  =  Z^  (r^)  =  Z^  (r^)  =  Z^  (r^)  =  0.  The  final  condition 

useful  for  checking  optimality  comes  from  the  Hamiltonian  value  condition,  which 
as  shown  in  (1 .16)  is  a  constant  negative  one  for  a  minimum  time  problem. 
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B.  MANEUVER  CHECKOUT  METHODOLOGY 


In  keeping  with  good  engineering  practice,  several  methods  of  V&V  are 
used  for  the  checkout  for  the  optimal  control  solutions.  The  first  of  these  is  to 
verify  that  the  magnitude  of  the  quaternion  vector  is  always  unity.  This  condition 
is  expressly  imposed  in  the  coding  of  the  problem  for  solving  in  DIDO.  Any 
deviation  of  the  magnitude  from  one  beyond  that  expected  from  machine 
rounding  would  be  an  immediate  red  flag  that  the  solution  is  invalid. 

The  second  method  of  V&V  is  an  inspection  of  the  Hamiltonian.  Since 
(4.1)  is  a  minimum  time  problem,  the  Hamiltonian  value  condition  and 
Hamiltonian  evolution  condition  specify  that  the  Hamiltonian  should  be  a  constant 
at  negative  one.  This  is  a  very  simple  and  powerful  check. 

The  third  V&V  method  is  to  verify  the  conditions  obtained  via  the 
transversality  equations.  The  transversality  equations  specified  that  the 
momentum  costates  are  equal  to  zero  at  the  conclusion  of  the  maneuver.  This  is 
another  necessary  condition  that  is  easy  to  check. 

The  final  and  most  involved  V&V  check  is  a  propagation  of  the  solution 
controls  using  a  simple  Simulink  model  of  the  LRO.  The  Simulink  model  will  use 
interpolated  control  signals  as  control  input.  The  LRO  dynamics  will  be  simulated, 
and  rates  integrated  over  time  to  produce  propagated  solutions  to  attitude, 
angular  velocity,  and  RW  momentum.  A  pictorial  view  of  the  Simulink  model  is 
shown  in  Figure  6.  Quaternions,  S/C  angular  velocities,  and  RW  momenta  with 
respect  to  time  will  be  plotted  and  checked  against  the  expected  results. 

In  addition  to  verifying  that  the  occupation  constraints  are  satisfied  for 
each  solution,  the  checkout  exercises  described  above  can  be  used  to  automate 
the  maneuver  checkout  process.  Suitable  tolerances  for  each  V&V  step  can  be 
imposed  by  the  S/C  operators.  If  these  tolerances  are  met,  then  the  maneuver 
can  be  considered  as  flight  ready  as  far  as  LRO  dynamic  feasibility  is  concerned. 
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Figure  6.  Simulink  propagator  used  for  maneuver  checkout. 

C.  INITIAL  RESULTS 

The  initial  optimal  control  problem  was  solved  using  DIDO  software.  The 
resulting  trajectories  of  the  primary  axes  of  the  S/C  are  shown  in  Figure  7.  The 
angular  distance  of  the  S/C  +z  axis  from  the  sun  is  shown  in  Figure  8.  The 
angular  velocity  of  the  S/C,  momentum  in  each  RW,  and  the  solved  control 
signals  are  shown  in  Figures  9,  10,  and  11,  respectively.  The  total  maneuver  time 
is  1,135  seconds,  which  represents  a  43  percent  improvement  over  the  status 
quo  approach  (see  Table  3). 
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Change  in  S/C  orientation  with  time 


1 


\ 


Inertial  X 

Figure  7.  S/C  orientation  evolution  for  initial  problem  formulation. 


0  200  400  600  800  1000 

time  (s) 

Figure  8.  Angle  between  S/C  +z  axis  and  sun  vector  for  initial  problem 

formulation. 
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S/C  Angular  Velocity 


Figure  9.  S/C  angular  velocity  evolution  for  initial  problem  formulation. 


RW  momenta 


Figure  10.  RW  momentum  storage  for  initial  problem  formulation. 
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Control  (h  dot) 


0  200  400  600  800  1000 

time  (s) 

Figure  1 1 .  Control  signals  for  initial  problem  formulation. 


As  seen  in  Figures  7  and  8,  the  optimal  control  problem  allows  the  S/C  +z 
axis  to  come  up  to  the  limits  of  a  pointing  violation,  without  violating  the 
constraint.  Comparing  with  the  direct  slew  and  dogleg  maneuvers  shown  in 
Figure  3,  the  improvement  appears  obvious.  While  the  solution  is  an  optimal 
solution  within  the  limits  of  computational  fidelity,  it  is  far  from  a  good  solution 
from  the  perspective  of  implementation.  The  control  signal  (see  Figure  11)  rapidly 
oscillates  between  the  maximum  positive  and  negative  values,  which  is  a 
possible  indication  that  active  control  is  not  needed  during  those  times  and  can 
be  zeroed  out.  This  in  turn  would  have  the  effect  of  smoothing  the  RW 
momentum  and  S/C  rate  trajectories. 

The  results  are  a  good  illustration  on  how  an  improperly  specified  problem 
can  lead  to  subpar  results.  Pseudospectral  methods  are  a  very  strong  tool  to  use 
in  solving  control  problems,  but  the  tools  are  useless  if  not  used  properly. 

A  readily  apparent  problem  is  that  most  of  the  changes  in  the  controls 
occur  during  a  small  fraction  of  time  at  the  beginning  and  end  of  the  simulation. 
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This  is  due  to  the  relatively  large  angle  through  which  the  S/C  must  slew 
compared  to  the  small  angular  velocity  limit.  When  DIDO  performs  its  algorithms, 
it  breaks  up  the  time  space  into  a  number  of  smaller  nodes  defined  by  the  user. 
With  increasing  node  count,  the  time  spent  on  calculations  increases 
exponentially  and  quickly  becomes  impractical.  For  this  simulation,  the  node 
count  was  50,  and  the  time  it  took  to  solve  on  the  computer  was  over  nine  hours. 
Not  only  is  this  an  impractically  long  time,  but  the  node  count  used  is  actually 
insufficient  for  a  proper  solution. 

The  solution  is  point-wise  and  must  be  interpolated  between  nodes.  That 
is,  the  solution  is  actually  discrete  and  not  continuous.  When  only  a  small  fraction 
of  the  time  space  is  devoted  to  meaningful  control,  these  discrete  time  pieces 
lack  the  resolution  to  capture  appropriate  changes  needed  for  a  good  solution. 

Another  issue  with  the  abovementioned  solution  is  the  non-smooth  nature 
of  the  RW  momenta  curves.  This  points  to  the  controls  being  possibly 
unnecessarily  actuated  in  opposite  directions.  Due  to  the  bound  on  control 
signals,  the  problem  formulation  yields  switching  functions.  If  the  switching 
function  over  time  is  oscillating  about  zero,  the  switch  can  for  all  practical 
purposes  be  set  to  zero.  That  is,  the  control  signal  can  be  set  to  zero.  A  section 
of  the  plot  of  the  switching  function  associated  with  the  first  control  signal  is 
shown  in  Figure  12.  In  this  portion,  the  switching  function  is  seen  to  take  on  very 
small  positive  and  negative  values.  If  these  are  presumed  to  be  zero,  simplifying 
the  control  signal  into  a  bang-off-bang  type  profile,  ostensibly  the  results  will  be 
better  behaved.  This  behavior  can  also  be  corrected  with  a  better  problem 
scaling. 
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Switching  function  S1 


xIO'® 


Figure  12.  Switching  function  for  first  control  for  initial  problem  formulation. 

D.  MODIFICATIONS  TO  OPTIMAL  CONTROL  PROBLEM 

The  combination  of  the  size  of  the  discrete  time  steps  and  the  angular 
velocity  limit  of  the  S/C  causes  the  solution  to  look  undesirable,  as  discussed 
above.  The  S/C  accelerates  in  rotation  to  its  limit  in  a  very  short  time  relative  to 
the  solution  time  space  and  then  coasts  most  of  the  way  until  it  decelerates  in  a 
very  short  relative  time.  Thus,  some  modifications  to  the  problem  specifications 
are  considered  to  allow  a  more  effective  approach  to  the  optimal  control  problem. 
Specifically,  the  issues  of  problem  scaling  and  the  specification  of  angular 
velocity  limits  are  examined.  The  switching  function  will  be  maintained,  as  the 
removal  of  the  angular  velocity  specification  will  help  remedy  this  issue. 

1.  Scaling 

The  main  problem  with  the  existing  problem  formulation  is  that  it  is  poorly 
scaled.  As  seen  in  (4.2),  many  product  terms  exist  between  various  state  and 
costate  variables.  This  means  that  if  the  magnitudes  of  these  variables  vary 
greatly,  the  problem  can  potentially  become  very  sensitive  to  slight  perturbations, 
leading  to  a  poor  solution. 
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A  method  to  mitigate  the  issue  of  poor  scaling  is  to  compare  the  state 
variables  to  their  corresponding  costate  variables.  By  roughly  matching  the 
orders  of  magnitudes  between  the  pairs  of  variables,  the  problem  becomes  better 
scaled.  An  additional  benefit  of  better  scaling  is  that  the  calculation  time  for 
obtaining  a  solution  is  significantly  reduced.  Plots  of  the  rate  costates  and  the 
Hamiltonian  for  the  problem  prior  to  scaling  are  shown  in  Figures  13  and  14.  The 
Hamiltonian  for  a  minimum  time  problem  should  be  negative  one.  However,  as 
seen  in  Figure  14,  the  actual  Hamiltonian  takes  on  values  that  deviate  from 
negative  one  by  orders  of  magnitude  including  a  peak  to  almost  14,000.  The 
costates  compared  with  their  corresponding  states  in  Figure  9  make  apparent  the 
disparity  in  order  of  magnitude. 
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Figure  13.  Rate  costates  for  initial  problem  formulation. 
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Figure  14.  Hamiltonian  for  the  poorly  scaled  problem  formulation. 

To  resolve  this  problem,  scaling  is  introduced.  When  a  problem  is  scaled, 
the  state  variables  are  multiplied  by  some  factor.  The  idea  is  to  introduce  factors 
such  that  when  the  process  is  followed  to  calculate  and  minimize  the  Hamiltonian 
and  develop  the  adjoint  equation,  the  equations  are  affected  by  changes  in 
variables  evenly.  That  is,  scaling  is  intended  to  prevent  a  condition  in  which  a 
change  in  a  variable  is  marginalized  due  to  it  being  orders  of  magnitude  different 
than  other  variables. 

For  initial  problem  formulation,  let  the  variables  Q,  Q,and  H  be  scaling 
factors  for  the  quaternions,  angular  velocities,  and  RW  momenta,  respectively. 
Then,  the  scaled  variables  are 

=  [q/Q,ro/r^,h/H] 

=  [q,®,h] 

By  substituting  these  scaled  variables  into  the  optimal  control  problem  and 
using  them  to  meet  the  necessary  conditions,  the  magnitude  of  the  costates  can 
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be  manipulated.  Consider  the  final  term  of  (4.2),  involving  the  fourth  quaternion 
costate.  Substituting  (4.4)  into  this  term,  we  obtain 

^  =  ...+1  +  +...  (4.5) 

From  (4.5),  it  is  apparent  that  the  order  of  magnitude  of  the  fourth 
quaternion  costate  will  vary  inversely  with  several  scaling  factors  for  quaternions 
and  angular  velocities.  Since  (4.1)  is  a  minimum  time  problem,  the  Hamiltonian 
must  equal  negative  one  throughout  the  maneuver.  By  correctly  choosing  the 
scaling  factors,  the  problem  becomes  better  scaled  such  that  a  small  change  in 
the  state  variable  does  not  cause  a  large,  undesirable  change  in  the  costate 
variables. 

The  procedure  to  find  the  appropriate  scaling  factors  is  in  some  sense  a 
trial  and  error  process.  Equation  (4.5)  represents  a  relatively  simple  portion  of  the 
Hamiltonian  in  terms  of  the  relationships  of  the  variables  with  each  other.  With 
many  other  terms  present  each  affecting  one  another,  it  is  not  possible  by 
inspection  to  determine  how  changing  the  scaling  of  one  variable  will  affect  all  of 
the  costates  and  in  turn,  the  Hamiltonian.  The  plots  of  states  versus  costates 
provide  a  starting  point  to  identify  disparity  in  orders  of  magnitude,  but  the  full 
effects  aren’t  apparent  until  the  entire  process  is  performed,  which  in  this  case  is 
a  successful  DIDO  run.  The  scaling  factors  that  are  used  for  the  problem  are 
shown  in  Table  4. 


Table  4.  Scaling  factors  for  initial  problem  formulation. 


Variable 

Scaling  factor 

Quaternions 

0.2 

Angular  velocities 

0.003 

RW  momenta 

4.0 
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The  resulting  state-costate  pairs  for  angular  velocity,  as  well  as 
Hamiltonian  are  shown  in  Figures  15,  16,  and  17.  The  rate  state  and  costate 
pairs  are  on  the  same  order  of  magnitude,  which  allows  for  a  higher  precision 
solution.  The  Hamiltonian,  while  not  cleanly  a  constant  negative  one,  is  orders  of 
magnitude  closer  to  the  desired  value  than  in  the  unsealed  problem.  The  scaled 
simulation  completed  running  in  3,922  seconds,  which  is  roughly  a  nine-fold 
improvement  over  the  unsealed  case.  The  results  are  further  explored  in 
Chapter  V.  It  is  important  to  note  that  the  scaling  factors  are  unique  to  each 
specific  problem,  and  those  in  Table  4  are  specifically  for  the  conditions 
discussed  in  this  chapter. 


S/C  Angular  Velocity 


Figure  15.  Initial  problem  formulation:  S/C  angular  velocities  for  scaled  problem. 
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Rate  costates 


Figure  16.  Initial  problem  formulation:  rate  costates  for  scaled  problem. 


Hamiltonian 


Figure  17.  Initial  problem  formulation:  Hamiltonian  for  scaled  problem. 

2.  Angular  Velocity  Limit 

Imposing  a  limit  on  angular  velocity  also  causes  the  developed  approach 
to  the  optimal  control  problem  to  yield  subpar  results.  There  are  several  ways  to 
mitigate  this  problem.  The  first  is  to  increase  the  number  of  nodes.  This  would 
allow  more  discrete  time  steps  to  capture  the  transients,  resulting  in  a  higher 
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fidelity  solution.  As  discussed,  this  quickly  takes  an  already  impractical  time 
required  for  calculations  and  turns  it  unmanageable,  and  is  not  a  good  path 
forward.  A  second  approach  is  to  break  up  the  problem  into  three  segments: 
acceleration,  coasting,  and  deceleration.  Continuity  of  the  S/C  attitude  and 
velocity  as  it  enters  and  exits  the  coasting  phase  could  be  presumed  using 
pseudospectral  knots  [16]. 

Another  mitigation  option,  and  the  path  chosen  for  this  thesis,  is  to 
eliminate  the  limit  on  angular  velocity,  as  it  is  a  constraint  that  is  necessary  only 
for  conventional  slew  design.  Using  the  optimal  control  approach,  the  limit  on 
angular  velocity  can  be  removed  since  the  capacity  of  the  RW  assembly  should 
drive  the  behavior  of  the  system.  In  doing  so,  the  solution  time  can  be  decreased. 
The  angular  velocity  limit,  while  necessary  for  conventional  approaches,  is  not 
needed  for  optimal  control. 

Additionally,  a  rough  calculation  was  done  to  find  the  maximum  S/C 
angular  velocity  if  it  were  to  do  a  74  degree  slew  with  a  bang-bang  control  profile, 
followed  by  a  128  degree  slew.  Whether  or  not  the  resultant  momentum 
accumulation  would  violate  a  RW  hardware  limit  is  evaluated.  Assuming 
maximum  control  torque  and  zero  accumulated  RW  momentum,  equation  (2.2) 
can  be  used  to  calculate  the  following: 


=  =[0.16, 0.16, 0.16, 0.16f/Vm 

d)  =  -/  '  ((Z)  X  {Jco  +  Zh)  +  Zh) 


d)(0)  =  -/  '  {cOq  X  (/(Z)q  +  Zh)  +  Zh) 


-0.0356 

-0.000842 

0.000802 


deg/ 


d)(0)|  =  0.0356  deg/ 


Then,  this  value  is  used  to  extrapolate  the  time  it  would  take  to  maneuver 
over  the  specified  angle  via  a  simple  double  integrator  function.  This  assumes  no 
gyroscopic  or  damping  effects. 
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(4.6) 


Equation  (4.6)  represents  the  angle  rotated  through  time  t  given  a 
constant  angular  acceleration  a,  assuming  the  S/C  is  initially  at  rest.  Thus,  to 
calculate  the  time  to  cover  a  given  angle  9  using  a  bang-bang  control  profile, 
one  half  of  the  angle  is  used  to  calculate  a  time  for  the  acceleration,  and  then  the 
time  doubled  to  account  for  the  deceleration  [15].  The  maximum  S/C  angular 
velocity  and  RW  momentum  accumulation  can  be  found  by  substituting  the 
longest  acceleration  time  into 


CO  =t  CO 

max  max 

h  =t  T 

max  max  max 


The  fidelity  of  this  calculation  is  not  very  high  considering  it  does  not  take 
into  account  the  effects  of  momentum  accumulation,  nor  does  it  apply  to  any 
specific  maneuver  other  than  the  subtended  angle,  but  it  does  give  a  data  point 
to  work  with. 


The  resulting  maximum  S/C  angular  velocity  was  2.1  degrees/sec,  which 
is  not  an  intuitively  unreasonable  value.  The  RW  accumulated  no  more  than  16 
percent  of  its  operational  limit  in  momentum,  and  the  maneuver  time  is  211  sec. 
While  a  very  simplified  calculation,  it  does  approximate  the  total  amount  of  time  it 
would  take  to  do  a  two  legged  maneuver  of  appropriate  magnitude. 

As  a  result,  the  rest  of  this  thesis  will  assume  no  limit  on  S/C  angular  rate, 
both  as  an  academic  exercise  and  as  an  exploration  on  what  the  LRO  is  capable 
of  if  its  wheels  are  utilized  to  its  design  capabilities.  Results  will  be  compared 
against  each  other  and  with  the  simple  estimate  developed  in  this  section. 
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V.  SUN  AVOIDANCE  MANEUVERS  FOR  LRO 


In  order  to  explore  the  utility  of  optimal  control  solutions  for  LRO  maneuver 
design  automation,  a  variety  of  different  operational  scenarios  were  considered. 
Specifically,  the  following  conditions  were  explored: 

1 .  All  RW  initially  have  zero  stored  momentum,  such  as  the  state  of 
the  S/C  immediately  after  a  momentum  dump. 

2.  An  initial  momentum  bias  exists,  due  to  accumulated  environmental 
torques. 

3.  Initial  net  stored  momentum  in  the  body  frame  is  zero,  but  individual 
RW  have  stored  momentum  due  to  null  redistribution  of 
momentum. 

4.  Additional  pointing  constraints  introduced  to  accommodate 
occultation  avoidance  of  other  sensors  and  instruments 

Within  conditions  two  and  three,  varying  magnitudes  of  initial  RW 
momenta  are  examined.  For  condition  three,  specification  (or  lack  thereof)  of  final 
RW  momentum  values  is  also  examined.  Finally,  automated  maneuver  checkout 
was  performed  for  each  of  the  results  using  the  methods  described  in  Section 
IV.B.  The  conditions  described  above  will  hereafter  be  referred  to  as  case  1, 
case  2,  etc. 

A.  CASE  1 :  ZERO  INITIAL  MOMENTUM  IN  ALL  WHEELS 

The  first  condition  explored  is  that  of  the  original  problem  statement, 
equation  (4.1),  as  modified  per  the  discussion  in  Section  IV. D.  There  is  no  stored 
momentum  in  any  of  the  wheels,  and  thus  each  wheel  has  the  capacity  to  store 
momentum  equal  in  magnitude  to  its  operational  limit  in  either  direction.  RWs  are 
typically  run  at  some  nominal  bias  rate  to  avoid  friction  effects  which  may  corrupt 
speed  measurements.  However,  these  bias  rates  are  typically  small  enough  to 
be  considered  zero  for  practical  purposes. 

The  trajectories  of  each  of  the  S/C’s  main  axes  are  shown  in  Figure  18. 
The  angle  that  the  S/C  +z  axis  makes  with  the  center  of  the  avoidance  cone  is 
shown  in  Figure  19. 
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Change  in  S/C  orientation  with  time 


Angie  between  SIC  Z  and  Sun  Vector 


Figure  1 9.  Angle  between  S/C  +z  axis  and  sun  vector  for  case  1 . 


The  paths  traced  by  the  S/C  axes  are  smooth  curves,  indicating  a 
maneuver  that  is  a  candidate  for  implementation.  This  is  in  clear  contrast  with  the 
haphazard  appearance  of  the  trajectories  obtained  using  the  original  unsealed 
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problem  in  Figure  7.  From  Figure  19,  it  is  clear  that  the  maneuver  comes  as 
close  to  the  keep  out  zone  as  possible  without  actually  violating  the  constraint. 
The  control  signals,  S/C  angular  velocities,  and  RW  momenta  follow  in  Figures 
20,  21,  and  22. 


Control  (\dot{h}) 


Figure  20.  Control  signals  for  case  1 . 


RW  momenta 


Figure  21.  RW  momenta  for  case  1 . 
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S/C  Angular  Velocity 


Figure  22.  S/C  angular  rates  for  case  1 . 


The  control  signals  are  consistent  with  the  bang-bang  maneuver  expected 
to  be  associated  with  switching  functions.  The  S/C  velocity  now  reaches  a 
maximum  of  1.3  deg/s  about  the  y  axis,  which  is  13  times  that  of  the  operational 
limits  used  on  the  LRO.  However,  even  though  the  angular  rates  are  larger  than 
the  0.1  deg/s  soft  limit,  the  RW  only  use  about  25  percent  of  their  operational 
momentum  capacity.  This  confirms  the  speculation  that  higher  rates  are  in  fact 
operationally  feasible. 

The  time  of  the  maneuver  is  189  seconds.  This  is  a  90  percent 
improvement  over  the  original,  rate  limited  dogleg  maneuver.  Compared  to  the 
estimate  made  in  Section  IV. D. 2,  it  represents  a  12  percent  improvement.  Thus, 
the  optimal  control  approach  succeeds  in  providing  a  method  for  automated 
maneuver  design,  while  at  the  same  time  improving  on  the  maneuver  time  over 
the  dogleg  approach. 

1.  Poor  Scaling  Example 

Another  interesting  solution  that  illustrates  the  need  for  proper  scaling  is 

presented  in  Figure  23.  In  this  case,  the  S/C  +z  axis  still  successfully  avoids  the 

obstacle,  as  required,  but  the  S/C  is  seen  to  rotate  in  the  opposite  direction.  The 
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resultant  maneuver  time  is  236  seconds,  which  is  25%  longer  than  case  1,  and 
11%  longer  than  the  estimated  dogleg  maneuver  time.  The  solution  is  a  feasible 
one,  and  optimal  within  the  numerical  fidelity  of  computation.  The  Hamiltonian  is 
minimized  for  this  given  problem  formulation,  number  of  nodes,  scaling,  etc.  This 
highlights  the  importance  of  determining  what  is  proper  scaling,  as  most 
indications  may  be  that  this  is  indeed  an  optimal  maneuver.  In  this  specific  case, 
there  was  a  disparity  in  the  magnitudes  of  state  costate  pairs.  This  results  in  a 
large  range  of  values  for  a  state  variable  which  could  have  little  effect  on  the 
Hamiltonian,  resulting  in  a  solution  that  is  only  locally  optimized. 


Change  in  S/C  orientation  with  time 


Figure  23.  S/C  attitude  evolution  for  poorly  scaled  case  1 . 

B.  CASE  2:  MOMENTUM  BIASED  CONDITION 

Case  2  explores  a  condition  in  which  all  RW  are  operating  at  the  same 
speed,  in  the  same  direction.  This  causes  a  net  stored  momentum  in  the  S/C 
body  frame,  called  a  momentum  bias.  Additionally,  the  initial  stored  momentum 
changes  the  maneuver  space  in  that  the  RWs  are  closer  to  one  limit,  such  as  the 
maximum  positive,  than  the  other.  If  the  fraction  of  the  operational  momentum 
limit  of  the  wheel  is  represented  by  r,  the  initial  conditions  for  the  wheels  is 
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written  as  (5.1).  Since  RWs  can  spin  in  either  direction,  r  can  take  on  values 
ranging  from  negative  one  to  positive  one. 


(5.1) 

Equation  (5.1)  is  used  as  part  of  the  specifications  of  the  optimal  control 
problem,  which  is  then  solved.  The  values  of  r  for  which  a  solution  was 
attempted  is  summarized  in  Table  5.  Table  5. 


Table  5.  Summary  of  initial  momentum  storage  values  for  case  2. 


Case  number 

r 

2.1 

+0.1 

2.2 

+0.25 

2.3 

+0.5 

2.4 

+0.7 

2.5 

-0.1 

2.6 

-0.25 

All  of  the  cases  contained  in  Table  5  are  solved,  and  referred  hereafter  as 
cases  2.1  through  2.6.  In  the  interest  of  brevity,  not  all  results  will  be  discussed, 
but  features  that  distinguish  the  results  will  be  summarized.  Most  of  the  cases 
essentially  build  on  each  other  as  a  confirmation  that  optimal  control  can  be  used 
to  automate  the  maneuver  design  process. 

The  first  case  of  interest  is  2.b,  where  r  is  +0.25.  This  value  for  r  will  also 
be  examined  in  case  3  and  case  4  for  continuity.  The  paths  that  the  S/C’s  body 
axes  trace  out  are  shown  in  Figure  24.  The  angular  distance  of  the  +z  axis  from 
the  center  of  the  obstacle  is  shown  in  Figure  25. 
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Change  in  S/C  orientation  with  time 


Figure  24.  S/C  attitude  evolution  for  case  2.2. 


Angle  between  S/C  Z  and  Sun  Vector 


Figure  25.  Angle  between  S/C  +z  axis  and  sun  vector  for  case  2.2. 

The  fact  that  a  momentum  bias  exists  in  the  system  tends  to  rotate  the 
entire  system,  and  this  is  evidenced  by  the  spiraling  motion  of  the  +x  and  +y 
axes.  This  corresponds  roughly  to  a  slow  rotation  about  the  +z  axis.  The  effect  is 
due  to  the  fact  that  as  initial  momentum  bias  is  increased,  the  a^iJw+Zh)  term 
in  equation  (2.2)  becomes  a  larger  fraction  of  the  total  torque  acting  on  the  S/C. 
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This  effect  becomes  larger  as  angular  velocity  and  initial  momentum  bias 
increase. 

Ostensibly,  there  is  some  value  of  r  at  which  the  problem  becomes 
infeasible  to  solve  due  to  the  cross  product  term  dominating  the  torque,  rendering 
control  signals  useless  at  small  angular  velocities.  If  dominating  is  defined  to 
mean  that  the  cross  product  term  is  10  times  the  magnitude  of  the  control  torque 
when  S/C  angular  velocity  is  0.5  deg/s  about  each  axis,  then  the  value  of  r 
which  would  cause  this  is  only  roughly  +/-0.17,  as  seen  in  Figure  26. 


Figure  26.  Ratio  of  Q)x{Jo)+Zr\^J  to  as  a  function  of  r,  with 

ft;  =  [0.5,0.5,0.5f . 

The  S/C  angular  velocities  and  momentum  storage  of  each  RW  are  shown 
in  Figures  27  and  28.  Compared  with  case  1,  the  S/C  attains  slightly  higher 
velocities.  The  RW  momenta  all  tend  in  the  same  direction,  spinning  down.  This 
behavior  is  responsible  for  the  spinning  motion  of  the  entire  S/C.  The  solution 
maneuver  time  is  224  seconds,  which  is  18  percent  longer  than  case  1,  and  six 
percent  longer  than  the  estimated  time  of  the  dogleg  maneuver. 
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S/C  Angular  Velocity 


Figure  27.  S/C  angular  rates  for  case  2.2. 


RW  momenta 


Figure  28.  RW  momenta  for  case  2.2. 


In  case  2.3,  where  r  is  0.50,  the  tendency  towards  spiraling  motion 
becomes  far  more  exaggerated,  as  seen  in  Figure  29.  Here,  the  S/C  appears  to 
spin  almost  out  of  control.  However,  all  boundary  and  path  conditions  are 
satisfied,  indicating  that  control  using  the  RWs  is  possible  despite  the  large 
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gyroscopic  torques.  Therefore,  an  operational  limit  on  the  allowable  momentum 
bias  should  be  specified. 


Change  in  S/C  orientation  with  time 


C.  CASE  3:  ZERO  NET  MOMENTUM  WITH  MOMENTUM  STORED  IN 
INDIVIDUAL  WHEELS 

Case  3  involves  a  zero  net  momentum  condition  as  with  case  1,  but  this 
time  the  RWs  will  have  some  accumulated  momentum  stored  in  them.  It  was 
determined  that  the  zero  momentum  condition  occurs  when  all  wheels  are 
operating  at  the  same  speed,  with  wheels  one  and  three  spinning  in  one 
direction,  and  wheels  two  and  four  in  the  other.  Using  the  same  notation  for 
fraction  of  maximum  as  with  (5.1),  the  initial  condition  of  the  wheels  can  be 
expressed  as 

(5.2) 
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A  summary  of  the  values  of  r  examined  for  case  3  is  shown  in  Table  6. 
Table  6.  Those  case  numbers  followed  by  an  “a”  denote  cases  in  which  a  final 
RW  momentum  condition  was  specified. 


Table  6.  Summary  of  initial  momentum  storage  values  for  case  3. 


Case  numbers 

r 

3.1  and  3.1a 

+0.1 

3.2  and  3.2a 

+0.25 

3.3 

+0.5 

3.4 

+1.0 

3.5  and  3.5a 

-0.1 

3.6  and  3.6a 

-0.25 

3.7 

-1.0 

1.  Case  3  with  Final  Momentum  Unspecified 

The  solution  for  S/C  attitude  motion  and  standoff  angle  of  the  Z-axis  from 
the  obstacle  are  shown  in  Figures  30  and  31 . 


53 


Change  in  S/C  orientation  with  time 


Inertial  X 

Figure  30.  S/C  attitude  evolution  for  case  3.2. 


Angle  between  S/C  Z  and  Sun  Vector 


Figure  31 .  Angle  between  S/C  +z  axis  and  sun  vector  for  case  3.2. 

The  most  important  observation  is  that  although  the  RWs  have  stored 
momentum  equal  in  magnitude  to  case  2.2,  the  fact  there  is  no  net  momentum 
bias  allows  the  solution  to  evolve  similarly  to  that  of  case  1.  In  fact,  the 
differences  in  the  path  between  case  1  and  case  2.2  are  so  minute  that  they  are 
not  discernable  from  the  figures.  The  S/C  angular  velocities  and  RW  momenta, 
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shown  in  Figures  32  and  33,  tell  a  similar  story.  Note  that  the  individual  RW 
momenta  in  Figure  32  trace  the  same  paths  as  in  Figure  21,  but  with  a  bias. 


S/C  Angular  Velocity 
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time  (s) 


Figure  32.  S/C  angular  rates  for  case  3.2. 


RW  momenta 


Figure  33.  RW  momenta  for  case  3.2. 


The  solution  maneuver  time  was  189  seconds,  which  is  the  same  as 
case  1 .  Cases  3.1 ,  3.5,  and  3.6  all  have  solutions  very  similar  to  that  of  case  3.2. 
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When  the  magnitude  of  initial  momentum  storage  in  individual  wheels  is 
increased  to  near  the  wheel  capacity,  the  results  change  dramatically.  Take  the 
extreme  example  of  case  3.7,  where  the  RWs  are  at  their  operational  limits  in 
terms  of  momentum  storage,  but  the  net  momentum  in  the  body  frame  is  still 
zero.  This  initial  state  of  the  RWs  severely  restricts  the  maneuver  space  of  the 
S/C  as  the  individual  wheels  cannot  spin  any  faster.  The  path  traced  by  the  body 
axes  and  the  angular  distance  of  the  +z  axis  from  the  sun  are  shown  in  Figures 
34  and  35. 


Change  in  S/C  orientation  with  time 
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Figure  35.  Angle  between  S/C  +z  axis  and  sun  vector  for  case  3.7. 

The  solution  time  is  283  seconds,  which  is  50  percent  greater  than  case  1 
and  37  percent  greater  than  the  estimated  time  for  the  dogleg  maneuver.  The 
S/C  +z  axis  traces  a  path  that  goes  around  the  obstacle  the  long  way  similar  to 
Figure  23,  but  this  time  due  to  the  restrictions  placed  by  RW  momentum  storage. 
The  RW  momenta  are  shown  in  Figure  36.  Note  that  unlike  in  Figure  33,  not  all 
wheels  are  being  actuated  simultaneously  at  all  times  due  to  the  operational 
limitations  placed  upon  them.  Case  3.7  proves  that  a  feasible  optimal  control 
maneuver  could  be  automatically  designed  even  in  a  very  challenging  situation. 
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RW  momenta 


Figure  36.  RW  momenta  for  case  3.7. 


2.  Case  3  with  Final  Momentum  Specified 

Cases  3.1,  3.2,  3.5,  and  3.6  were  also  solved  while  enforcing  a  final 
boundary  condition  on  the  RW  momenta.  These  cases  will  be  referred  to  as 
cases  3.1a,  3.2a,  3.5a,  and  3.6a.  At  times,  it  may  be  necessary  to  drive  the  RWs 
to  a  desired  state,  such  as  to  redistribute  stored  momentum. 

To  calculate  the  final  momentum  specification,  the  fact  that  system 
momentum  is  conserved  must  be  used.  Since  the  S/C  is  stationary  at  the  start 
and  the  finish  of  the  maneuver,  it  follows  that  all  of  the  system’s  stored 
momentum  must  be  entirely  in  the  RW.  The  momentum  stored  in  the  RW 
expressed  in  x,  y,  and  z  components  can  be  found  by  left  multiplying  the  RW 
momentum  vector  by  the  Z  matrix.  This  value  must  be  the  same  at  the 
beginning  and  end  of  the  maneuver.  Since  the  S/C  orientation  changes,  the  Z 
matrix  in  the  inertial  frame  is  not  the  same  at  the  beginning  and  end  of  the 
maneuver,  and  thus  must  be  transformed  using  a  DCM.  Let  C  be  the  DCM 
described  by  (2.4),  corresponding  to  a  transformation  from  =  [0,0,0, if  to  a 
desired  quaternion.  Then,  the  desired  boundary  condition  can  be  found  as 
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CfZh^  =  C^Zh^ 


(5.3) 


h^-  =  Z^Cj^CgZh^ 

In  (5.3),  the  dagger  denotes  the  pseudoinverse  of  the  Z  matrix.  Since  Z 
is  not  a  square  matrix,  it  is  not  invertible.  However,  a  pseudoinverse  can  be  used 
to  solve  the  equation.  The  pseudoinverse  of  Z  is  calculated  using  equation  (5.4) 
where  Z*is  the  Hermitian  conjugate  of  Z  [17]. 

Z^=(Z*Z)‘Z*  (5.4) 

Equation  (5.3)  is  imposed  as  the  final  boundary  condition  on  RW 
momentum.  This  boundary  condition  also  serves  as  a  proof  of  concept  in  the 
ability  to  despin  wheels  while  performing  a  maneuver.  In  this  way,  momentum 
redistribution  can  be  performed  during  a  reorientation  maneuver.  The  S/C 
orientation  evolution  and  angle  between  +z  axis  and  the  obstacle  are  shown  in 
Figures  37  and  38. 


Change  in  S/C  orientation  with  time 


Figure  37.  S/C  orientation  evolution  for  case  3.2a. 
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Figure  38.  Angle  between  S/C  +z  axis  and  sun  vector  for  case  3.2a. 

The  paths  traced  out  by  S/C  axes  in  Figure  37  are  slightly  but  noticeably 
different  than  those  in  Figures  30  and  31.  The  maneuver  takes  210  seconds, 
which  is  11  percent  longer  than  case  1  and  case  3.2,  and  approximately  equal  in 
length  to  the  estimated  time  for  a  dogleg  maneuver.  By  forcing  the  RW  to  a 
desired  end  state,  an  additional  cost  is  incurred  in  the  form  of  a  longer  maneuver, 
but  there  may  be  advantages  operationally  for  increasing  this  extra  cost.  The  RW 
momenta  are  shown  in  Figure  39.  It  can  be  verified  that  the  final  condition  of  the 
RWs  are  that  they  are  not  spinning,  the  condition  which  is  consistent  with  the 
application  of  (5.3). 
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RW  momenta 


Figure  39.  RW  momenta  for  case  3.2a. 


D.  CASE  4:  MANEUVERS  WITH  ADDITIONAL  OCCULTATION 

AVOIDANCE  CONSTRAINTS 

Cases  1,  2,  and  3  explored  a  maneuver  in  which  there  was  only  one 
pointing  constraint  which  needs  to  be  satisfied.  In  reality,  as  discussed  in  Table 
1,  there  are  multiple  constraints  for  the  LRO.  To  ensure  that  a  maneuver  with 
multiple  occultation  constraints  could  be  designed,  two  additional  constraints  are 
introduced  to  the  problem  formulation.  The  first  additional  constraint  is  a  26.5 
degree  cone  simulating  the  Earth’s  disc  plus  a  15  degree  standoff  range  that  the 
S/C  +x  axis  (simulating  the  bore  sight  of  a  star  tracker)  must  avoid.  The  second 
obstacle  is  a  25.5  degree  cone,  centered  at  the  sun  and  therefore  oriented  in  the 
same  direction  as  the  63  degree  cone  from  previous  cases.  For  this  obstacle,  the 
S/C  +y  axis  (another  simulated  star  tracker  bore  sight)  must  avoid  the  cone.  The 
two  new  pointing  constraints  that  must  be  added  to  the  optimal  control  problem 
formulation  are: 


6^  =  cos  ^{E»x)  >  26.5 
6^.  =  cos  '(S’y)  >  25.5 


(5.5) 
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In  (5.5),  E  is  the  unit  vector  pointing  in  the  direction  of  the  Earth,  and  x 
and  y  are  unit  vectors  pointing  along  the  S/C  +x  and  +y  axes,  respectively. 
These  obstacles  are  placed  such  that  both  constraints  will  be  violated  by  the 
paths  traced  in  case  3.2,  as  shown  in  Figure  40. 


Change  in  S/C  orientation  with  time 


Figure  40. 


Inertial  X 

Constraint  violations  using  the  path  solved  for  case  3.2. 


In  Figure  40,  the  additional  keep  out  cones  are  colored  such  that  they 
correspond  to  the  colors  of  the  axis  which  must  avoid  the  particular  cone.  The 
S/C  +x  axis  (blue)  clearly  violates  the  blue  cone,  and  the  S/C  +y  axis  (green) 
clearly  violates  the  green  cone.  The  standoff  angles  of  these  axes  to  their 
respective  obstacles  are  shown  in  Figures  41  and  42,  verifying  the  pointing 
violation. 
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Figure  41 .  Angle  between  S/C  +x  axis  and  Earth  vector  for  case  3.2. 


Angle  between  S/C  Y  and  Sun  Vector 


Figure  42.  Angle  between  S/C  +y  axis  and  sun  vector  for  case  3.2. 

A  new  maneuver  solution  was  obtained  after  enforcing  the  additional 
occultation  constraints.  The  overall  S/C  orientation  and  the  angular  distances  of 
each  axis  from  their  respective  obstacles  are  shown  in  Figures  43,  44,  45,  and 
46. 
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Inertial 


Change  in  S/C  orientation  with  time 


Figure  43.  S/C  orientation  evolution  for  case  4 


Angle  between  S/C  Z  and  Sun  Vector 


Figure  44.  Angle  between  S/C  +z  axis  and  sun  vector  for  case  4 
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Angle  between  S/C  X  and  Earth  Vector 


Figure  45.  Angle  between  S/C  +x  axis  and  Earth  vector  for  case  4 
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Figure  46.  Angle  between  S/C  +y  axis  and  sun  vector  for  case  4 

The  above  figures  demonstrate  that  automated  maneuver  design  is 
feasible  with  the  introduction  of  additional  occultation  constraints.  The  most 
remarkable  result  is  that  the  solution  time  is  192  seconds,  less  than  two  percent 
longer  than  case  3.2  with  only  one  obstacle.  It  is  unlikely  that  this  level  of 
performance  would  be  obtained  by  manually  designing  a  slew  as  per  the 
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approach  outlined  in  Chapter  III.  The  S/C  angular  velocities  and  RW  momenta 
are  shown  in  Figures  47  and  48,  illustrating  the  feasibility  of  the  maneuver. 


S/C  Angular  Velocity 


Figure  47.  S/C  angular  velocities  for  case  4. 


RW  momenta 


Figure  48.  RW  momenta  for  case  4. 
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E. 


MANEUVER  CHECKOUT 


Maneuver  checkout  is  done  in  accordance  with  the  methods  described  in 
section  IV. B.  While  the  checkout  steps  were  performed  for  all  cases,  the  results 
for  case  4,  the  most  complex  case,  are  presented  here. 

The  first  check  is  to  verify  that  the  quaternion  norm  condition  is  satisfied 
throughout  the  maneuver.  The  quaternion  norm  error  with  respect  to  time  is 
shown  in  Figure  49.  This  condition  is  clearly  satisfied. 


^.|q-4  Quaternion  norm  error 


Figure  49.  Checkout  of  the  quaternion  norm  for  case  4. 


Next,  the  Flamiltonian  is  evaluated.  For  a  minimum  time  maneuver,  the 
Flamiltonian  should  be  a  constant  at  negative  one.  The  Flamiltonian  for  case  4  is 
shown  in  Figure  50.  While  at  first  it  may  appear  that  the  Hamiltonian  varies  wildly, 
note  the  scale.  Tiny  variations  such  as  those  seen  in  Figure  50  are  expected  and 
are  due  to  numerics.  The  variation  of  the  Hamiltonian  is  4.78x10'^.  Therefore,  it 
can  be  concluded  that  the  Hamiltonian  value  condition  is  satisfied. 
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Hamiltonian 


Figure  50.  Checkout  of  the  Hamiltonian  for  case  4. 


Additionally,  the  boundary  condition  imposed  by  the  transversality 
condition  is  verified.  For  case  4,  the  momentum  costates  at  the  end  of  the 
maneuver  should  be  equal  to  zero.  The  actual  momentum  costate  functions  are 
plotted  in  Figure  51.  While  not  exact  due  to  numerics,  the  momentum  costates  at 
the  end  of  the  maneuver  are  close  to  zero.  Thus,  the  transversality  conditions  are 
also  satisfied. 
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Momentum  costates 


Lastly,  the  control  signals  are  interpolated  and  propagated  using  a 
Simulink  propagator.  The  propagated  quaternions,  S/C  angular  velocities,  and 
RW  momenta  should  closely  approximate  those  obtained  from  DIDO.  The 
propagated  results  are  shown  in  Figures  52,  53,  and  54.  The  results  show  that 
the  solution  control  signals  result  in  the  expected  S/C  behavior.  Figure  55  shows 
the  propagated  avoidance  angle  for  all  three  axes. 
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Figure  52.  Attitude  checkout  for  case  4. 


S/C  Angular  Velocity 
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Figure  53.  Angular  rate  checkout  for  case  4. 
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RW  momentum  (N-m) 


RW  momentum 
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Figure  54.  RW  momentum  checkout  for  case  4. 
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Figure  55.  Standoff  angles  checkout  for  case  4. 
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As  shown,  all  checkout  exercises  for  case  for  were  satisfactory.  This 
shows  that  the  maneuver  is  dynamically  feasible  and  can  be  implemented  in 
flight.  For  actual  implementation,  an  approach  similar  to  that  used  for  optimal 
control  implementation  on  the  TRACE  S/C  can  be  used  [14].  Note  also  that  other 
constraints,  such  as  thermal  constraints  can  be  added  as  part  of  the  optimal 
control  problem  formulation  so  that  the  design  and  checkout  process  can  be 
automated  for  these  as  well. 
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VI.  CONCLUSIONS  AND  FUTURE  WORK 


Optimal  control  theory  provides  the  tools  necessary  to  formulate  and  solve 
control  problems  to  optimize  control  paradigms.  The  LRO,  which  uses  a  manual 
approach  to  maneuver  design,  stands  to  benefit  from  automation  of  maneuver 
design.  To  meet  the  thesis  objectives,  a  slew  maneuver  was  designed  such  that 
without  intervention,  a  pointing  constraint  violation  would  occur.  The  LRO  would 
normally  use  a  dogleg  maneuver  to  prevent  this  pointing  constraint  violation. 
When  the  problem  was  framed,  formulated,  and  solved  as  an  optimal  control 
problem,  reduction  in  maneuver  time  of  up  to  12  percent  was  realized.  If  the 
conventional,  rate  limited  maneuvers  are  considered,  the  time  savings  is  up  to  90 
percent.  Additionally,  the  optimal  control  approach  was  able  to  design  maneuvers 
in  challenging  conditions,  such  as  an  initial  momentum  bias  in  the  system  or  the 
introduction  of  additional  pointing  constraints,  without  human  intervention. 
Optimal  control  thus  provides  an  automated  approach  to  the  design  and 
checkout  of  maneuvers,  simplifying  ground  segment  operations  and  increasing 
the  potential  for  science  data  collection. 

Future  work  can  include  several  topics.  One  topic  is  the  extension  towards 
automation  and  prioritization  of  consecutive  maneuvers.  It  may  be  that  when 
multiple  maneuvers  are  desired,  such  as  multiple  science  collects,  that  there  is 
benefit  in  considering  a  specific  order  for  them.  Another  topic  for  future  work  can 
be  to  develop  an  approach  to  determine  an  optimal  paradigm  for  momentum 
dumping,  similar  to  the  one  demonstrated  on  the  ISS  [18].  Such  an  approach 
could  minimize  the  number  of  times  momentum  dumping  occurs,  or  combine 
momentum  dumping  with  other  slew  maneuvers  to  increase  the  time  utilized  for 
science  data  collection  and/or  save  on  fuel  resources. 
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